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This paper discusses optimization tech- 
niques for the determination of com- 
plex permittivity and permeability in 
transmission lines. The traditional theo- 
retical model using scattering parame- 
ters is extended into a mathematical 
regression model that can be solved 
with widely accepted numerical tech- 
niques. This new model produces accu- 
rate primary mode results for the 
samples tested including nonmagnetic 
and magnetic materials with high 
dielectric constants. An extension of 
the model includes responses due to 
higher order modes. The general model 
determines parameters to specify the 



spectral functional form of complex 
permittivity and permeability and is ca- 
pable of small corrections to indepen- 
dent variable data including angular 
frequency, sample length, sample posi- 
tion, and cutoff wavelength. TTie 
method provides reliable determination 
for both low and high permittivity ma- 
terials. 
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1. Introduction 

A constrained nonlinear optimization procedure 
is presented for the determination of complex per- 
mittivity and permeability spectra from scattering 
parameter (5 -parameter) data taken from an auto- 
matic network analyzer (ANA). The procedure has 
been used successfully for reliable characterization 
of permittivity and permeability of many different 
test samples. In addition, it provides a basis for the 
analysis of multi-mode field data and for the deter- 
mination of experimental systematic uncertainty. 

Previous work in this area involved the determi- 
nation of permittivity and permeability on a point- 
by-point basis with explicit or implicit solution of a 
system of nonlinear scattering equations at each 
particular frequency (see [1, 5]). Inaccurate results, 
however, may arise when numerical singularities 
occur at frequencies corresponding to integer multi- 
ples of one half wavelength of the material. The 



same system of nonlinear scattering equations is 
used in this study. Here they are solved in the sense 
of least squares over the entire range of measure- 
ment by determining the best Laurent series ap- 
proximations to permittivity and permeability 
consistent with linearity and causality constraints. 
Points of singularity may be de-emphasized to 
lessen the effect of highly uncertain data points. 
This effort determines complex permittivity and 
permeability from two-port 5 -parameter data using 
the primary (or fundamental) mode field behavior 
in various materials. Physical measurements of the 
5 -parameter data are made with an automatic 
network analyzer. Fundamental mode 5 -parameter 
relationships are then used to solve for permeabil- 
ity and permittivity, ^i-{(l}) and e(w) respectively, as 
a function of the angular frequency <u. 
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The general functional form for fj, and e was de- 
cided upon after evaluating several different poly- 
nomial and trigonometric relations. The evaluation 
criteria for the functional form for permittivity and 
permeability was based on the reproducibility of 
the 5 -parameter data in terms of the minimum to- 
tal least square approximation to these data. The 
best overall approximation was selected. The gen- 
eral functional form for /x and e with acceptable 
results involves only the first two terms of a Lau- 
rent series. 

The optimization approach in this research is an 
implicit function regression model. This model is 
solved by the orthogonal distance regression pack- 
age ODRPACK [2]. This estimation package al- 
lows for adjustments in input parameters to 
compensate for measurement uncertainties. Ad- 
justments here are limited to the sample length, 
sample position in the waveguide measurement fix- 
ture, and the cutoff wavelength. Orthogonal dis- 
tance regression is intended to compensate for 
slight uncertainties in the independent variable 
(angular frequency, co) as well as the dependent 
variable (observed 5-parameter data). 

The physical model is outlined in Sec. 2, with the 
related mathematical model discussed in Sec. 3. 
Numerical considerations are covered in Sec. 4 and 
conclusions and future directions discussed in Sec. 
5. 

2. Scattering Parameter Relations 

The equations described below relate the mea- 
sured two-port scattering parameters (5 -parame- 
ters) to the permittivity and permeability of the 
material. First, in order to develop the scattering 
equations, the following notation is used. Let 

^r-J£r \£o=£r£o 
and 

fJ. = [fJ-R -jfJ. r]iJ.O = IJ.RIJ.0, 

the permittivity and permeability of a sample mate- 
rial, where eo and /xo are the permittivity and per- 
meability of a vacuum, and e^ and /mr are the 
relative complex permittivity and permeability. 
Next let Cvac and ciab be the speed of light in a vac- 
uum and the laboratory, respectively, and for a 
given frequency/, let 

co = 2Trf, 



be the corresponding angular frequency. Then 



■ ^ vac ^ '*■ c I ' 



and 



70 



=W (^)'-©' 



represent the propagation constants in the material 
and air, respectively, where 7 =-\/—l and A^ is the 
cutoff wavelength in the waveguide measurement 
fixture, where the subscript 1 refers to the funda- 
mental mode. The expression for the transmission 
coefficient z is 

2;=exp(-yL), 

where L is the sample length. The reflection coeffi- 
cient is 



r= 



Mo y 

no y 



or 



Clab 



' £r 



Clab 



y 



* 

£r 



+ 1 



for coaxial line when -: > 0. 

Aci 



It is assumed that the total length of the sample 
holder is 



Laic — L + L1 + L2 



(1) 



as shown in Fig. 1, where Li and Li are the distances 
from the calibration reference planes to the sample 
faces for ports 1 and 2, respectively. 

For a two-port device the expressions for the 
measured 5 -parameters are obtained by the solu- 
tion of a related boundary value problem. The ex- 
plicit expressions for the scattering relations of the 
fundamental mode are assumed to functions of Ac, 
are given by 
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Figure 1. Dielectric sample in waveguide. 



Sn(\cd=Ri[^f^], (2) 

S,2(Kd=RiR2['-{zirf^]' (3) 

S2i(K)=R^R2[ \^^~rj ], (4) 

Sn(X.d=Rl[j^^l (5) 

where Ri = exp( — yo Li) and i?2= exp( — yo L2) are 
the reference plane transformation expressions. 

3. The Mathematical Model 

3.1 fitico), e*(io), and the Regression Model 

The mathematical problem can be stated as that 
of finding parameters to a prespecified functional 
form for the complex functions iir{o>) and e^(w). 
5 -parameter data acquired from the ANA for se- 
lected samples of materials provide complex values 
for each of 5ii, 52i, Sn, and ^22, at n different fre- 
quencies ranging from 1 to 18 GHz. The model de- 
termines the parameterization of /t| and e^ that 
best reproduces simultaneously the four 5 -parame- 
ters for the n observations, given the sample data, 
two reference plane positions, and the sample 
length. 

The general form for /A^(a)) and e|(w) uses 
terms from the Laurent series: 



f{o>)= 1 



tti 



{l + biOi)'' 



(6) 



where ai and bi are complex scalars in the /th term. 
The solution procedure here is an implicit function 
regression and corresponds to the parameterization 
of two terms in a truncated Laurent series. Note 
that Eq. (7) automatically satisfies Kramers-Kronig 
relations (see, for example, [1] or [4]) for disper- 
sion. The functional form for ij^r{o>) and e|(w) 
follows: 



/A^(ft>) =f{,(a) = 



eH<^)^f(<o) = 



Ai 



+- 



1+Bio) (l+B2(o) 



■ + 



l+BiO) (l+Bjcaf 



For notational convenience, the terms of the second 
truncated Laurent series ^1,^2, Bi, and B2 will be 
referred to as A3, A4, B3, and B4. 

The solution procedure to determine the com- 
plex parameters At , Bi is the minimization of the 
sum of the squares of the uncertainties between the 
predicted and observed 5 -parameters, 

min( i) I Si\-Pi\ i' + i 5i!-i'i*2i' + 

\l/2 
l52*l-i'2*i'+l5^-P22i'j (7) 



where, for ij = 1, 2, 5,* represents the fcth observed 
Sij scattering parameter at frequency m and 

Pij =Sij(Xci, (>>k, Ai, . . . , A4, Bi, . . . , -64), 

is the corresponding predicted scattering param- 
eter. We let \z\ = y^di(zy + ^(zf represent the ab- 
solute value of a complex scalar z. Note that Eq. (8) 
is equivalent to the minimization of the sum of the 
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squares of the real and imaginary parts of each of 
the 5-parameters, i.e., 



n 2 2 



min(i 2 2 ^(si-piy+^isti-piy) 



1/2 



-Jt=l ,=1 j=l 



(8) 

The model reported in Sec. 4.4 results uses the for- 
mulation in Eq. (8). 

3.2 Adjustments to the Model Inputs 

In standard ordinary least squares regression 
models, the observed responses are assumed to 
contain some uncertainties either produced by the 
phenomena under examination or introduced by 
the device that measures the events. In addition, 
certain independent and dependent variable pairs 
may not be as reliable as others due to an in- 
creased variance in the uncertainties in the depen- 
dent variable for particular values of the 
independent variable. The application of a reduced 
weight to points of questionable reliability provides 
the modeler with a means to de-emphasize such 
points to find a more appropriate regression solu- 
tion. 

Usually the independent variable can be con- 
trolled, and the precise value of each of the obser- 
vations is well known. An orthogonal distance 
regression model provides the modeler with the ad- 
ditional ability to assume that the independent 
variable, in this case frequency, may contain some 
uncertainty as well. Allowances for this type of un- 
certainty can, in some cases, greatly improve the 
approximation. For this particular model and the 
samples tested in this study, the uncertainty in the 
independent variables is sufficiently small to allow 
the modeler to assume that an ordinary least 
squares approximation provides an adequate solu- 
tion. 

Other model parameters such as sample length, 
sample position in the waveguide, and cutoff wave- 
length are sufficiently sensitive to require slight 
perturbation. Each of these inputs is required by 
the 5 -parameter equations and is considered to be 
known by the modeler. Since these inputs are not 
always known exactly, each may be perturbed 
slightly, as determined by the problem solver, to 
improve approximation. This allows the user to ad- 
just for measurement uncertainty. 

In particular, since incorrect specification of the 
sample position, Li, in waveguide affects the value 
of phase in the reflected 5 -parameter data, an ad- 
ditional parameter, Pn, is included in Ri, 



Similarly, uncertainties in Lair and L2 are repre- 
sented by analogous means with parameters jSmr 
and /3i2, respectively. With Eq. (1), the total length 
L of the sample is completely determined by 

L=(Lair+/3L.i,)-(Ll + /3L,+L2+/3L2) 

and is parameterized by the values of Pl^„ Pli, and 

Another parameterized correction, )3a, is in- 
cluded in the cutoff wavelength Ac , as Ac , + jSa . In- 
accuracies in the waveguide dimensions due to the 
milling process can affect the value for the cutoff 
wavelength. In addition, the higher order model 
discussed below requires an additional wavelength 
cutoff value, Xc2. for a higher mode solution. Since 
it is not known a priori which higher modes will be 
present, a variation in Ki from 0< Ac2< Aci + /3a, 
enables the solution procedure to find the value for 
Ac2 that best improves the 5 -parameter approxima- 
tion for the higher order terms. 

One can formulate an invariant model with re- 
spect to reference planes for the problem discussed 
here, as suggested in [1], to remove uncertainty in 
L 1 and in L2. This approach uses Eq. (3) or (4), and 
the determinant of the 5 -matrix, 

det(5) = SuSii - SuS2i, 

and solves for the parameters Ai,. . . , Aa, Bi . . . , 
54. Simplification of det(5) yields a formula with 
neither Li nor Li. This reduced model remains de- 
pendent on a precise knowledge of both Lair and L 
and the cutoff wavelength, A^. This approach can 
be used interchangeably with the original model 
that contains Li and L2. The original model uses 
twice the number of observations and seems to 
produce more accurate approximations to the 5- 
parameter data. 

3.3 Higher Order Modes 

In samples with a high dielectric constant, the 
observed 5 -parameter data may exhibit responses 
due to modes other than the fundamental mode. 
These responses are the result of resonances of the 
higher mode in the material. Earlier work does not 
include this information in the computation of /n« 
and er as higher modes are ignored. In this section 
we describe an enhancement to the previously de- 
fined model that does include higher order re- 
sponse data. 

Because of the similarity between the response 
of the primary mode and the higher modes (see 
Sec. 4.4), an additional term is added to each of the 
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four equations (see Eqs. 2-5). This numerical 
model includes higher order mode structure in the 
5 -parameter approximation and is a simple exten- 
sion of the primary mode model described in 
Sec. 2} The terms that approximate the higher 
order modes are identical to the primary mode 
term with the exception that each term is scaled. 
For the evaluation of the higher order terms, all 
input parameters are unchanged except for the 
parameter for the higher mode cutoff wavelength, 
Ac 2, that is allowed to decrease. 

The explicit forms of the new set of equations 
use the Eqs. (2-5), and are denoted as Dn, Dn, Diu 
and D22- Associated with each higher order term is 
a scaling parameter, ySi, . . . , yS4. In particular, to 
include the primary mode and one higher-order 
mode, the following approximate model is used: 

I>ii = 5ii(AcO + ySi5ii(Acj), (9) 

i)21=52l(A.,) -t- /3252l(A<:0, (10) 

I>12=5l2(A<:0 + yS35l2(Ac2), (11) 

I>22 = 522(Ac.) -f PaS22{\c^), (12) 

where Ac 2 is determined by ODRPACK. 

3.4 The Initial Solution Procedure 

In the Nicolson-Ross-Weir procedure (5, 7) the 
equations for the scattering parameters are com- 
bined to allow the system of equations to decouple. 
This decoupling yields an explicit equation for the 
permittivity and permeability as a function of the 
scattering parameters on a point-by-point basis. 
This solution procedure is the basis of the com- 
monly used techniques for obtaining permittivity 
and permeability. Unfortunately, these equations 
are not well-behaved for low-loss materials at fre- 
quencies that correspond to integer multiples of 
one half wavelength in the sample. 

The Nicolson-Ross-Weir procedure that is im- 
plemented provides a good initial approximation to 
A\, A3, B\ and Bi. All other parameters are initial- 
ized to zero. The estimated values for permittivity 
and permeability are determined on a point-by- 
point basis by frequency. The corresponding scat- 
tering parameters are computed with these values 
and then compared to the observed values. The 
computed values for \it and cr that provide the 
closest agreement between the observed and pre- 
dicted 5 -parameter data are used as the initial val- 
ues to the regression model. 



' The justification of new terms added to each of Eqs. 2-5 is 
based solely on empirieal evidenee found when solving for the 
primary mode. 



4. Numerical Considerations 

4.1 ODRPACK Orthogonal Distance Regression 
Package 

Briefly, ODRPACK [2] is an implementation of a 
trust region Levenberg-Marquardt algorithm. This 
type of trust region approach adaptively determines 
the region in which the linear approximation closely 
resembles the nonlinear model. The procedure al- 
lows both an ordinary least squares model, in which 
the uncertainties are assumed to be only in the de- 
pendent variable, and, an orthogonal distance re- 
gression model, where uncertainties also exist in the 
independant variables. 

First order derivatives for the Jacobian matrices 
can be numerically approximated (finite difference 
approximation), or can be user-supplied analytical 
derivatives. The procedure performs automatic 
scaling of the variables if necessary, as well as deter- 
mination the precision of the model in terms of ma- 
chine precision. The ODRPACK includes many 
other features that assist the user in the modeling 
process. The model automatically determines the 
number of digits in the model, checks analytical 
derivatives provided by the user, and automatically 
selects many of the input parameters for the user, if 
desired. 

Iterations are stopped in ODRPACK when any 
one of three stopping criteria is met. Two of these 
indicate that the iterations have converged to a so- 
lution. Sum-of-squares convergence indicates that 
the change in sum-of-squares observational uncer- 
tainty is sufficiently small. Parameter convergence 
indicates that the change in the estimated parame- 
ters is sufficiently small. The third stopping crite- 
rion is a limit on the number of iterations. 

4.2 Initial Conditions 

Many of the input options for ORDPACK can be 
set to their default values, as is done in this model. 
The most significant input parameters for modeling 
permittivity and permeability are the initial values 
for Ai, A3, Bi and B3. Sensitivity to the initial solu- 
tion for these parameters is discussed below and the 
selection of initial settings is covered in Sec. 3.4. 
When higher modes are included, the solution for 
the primary mode is used as the initial guess for the 
higher mode model. For all of the parameters that 
define the physical model except those for ^lR and 
er mentioned above, the standard laboratory values 
are used. All additional parameters are initialized 
to zero. 
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4.3 Sample Characteristics 

A total of seven samples were modeled to deter- 
mine permittivity and permeability with ANA two- 
port 5 -parameter data. The sample characteristics 



appear in Table 1. Figure 2 contains selected plots 
of the observed ANA 5 -parameter data versus 
frequency for sample 6. 



Table 1. Sample characteristics 





Material 




Length 


Li 


L2 


Initial 


Initial 




name 


Identifier 


(m) 


(m) 


(m) 


e^ 


M^ 


1 


Cross-linked polystyrene 


rexa240889 


2.407x10-2 


0.000 


1.3208x10-3 


(2.53,0.002) 


(1.0) 


2 


1723 glass 


172a240889 


1.015x10-2 


0.000 


1.5250x10-2 


(6.15,0.04) 


(1.0) 


3 


Loaded polymer 


112a050290 


2.540x10-2 


0.000 


0.0000 


(5.75,0.23) 


(1.6,0.1) 


4 


Yttrium iron garnet 


YIG 


1.766x10-2 


0.000 


7.7380x10-3 


(10,0.2) 


(1.0.4) 


5 


Nickel ferrite 


ttlal20490 


1.013x10-2 


0.000 


1.5263x10-2 


(11.5,0.1) 


(0.85,0.1) 


6 


Barium titanate mix 1 


barium 


7.632x10-3 


0.000 


1.7768x10-2 


(265,1) 


(1.0) 


7 


Barium titanate mix 2 


fircnch 


2.427x10-2 


0.000 


1.1252x10-3 


(105,1) 


(1.0) 
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Figure 2. Barium titanate mix 1 (a) 9i [5ii], (b) 5 [5ii], (c) SR [52i], (d) 3 [^21]. 
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4.4 Numerical Results 

The first results reported are those for the pri- 
mary mode model for all seven samples described 
in Sec. 4.3. The first set of plots (Fig. 3) include the 
observed 5 -parameter data (dots) for the YIG 
sample (sample 4) from the ANA overlaid by the 
predicted data (line) found by the model. The cor- 
responding residual plots^ for this sample appear in 
Fig. 4. 

For the first four samples the predicted and ob- 
served data are nearly identical. The residual plot 
for cross-linked polystyrene (sample 1) as shown in 
Fig. 5 reveals the systematic uncertainty due to the 
ANA. 

Fig. 6 illustrates a 5 -parameter primary mode 
solution for barium titanate mix 1 with a high 
dielectric constant and its corresponding residual 



plot. For three of the samples with high dielectric 
constants, the new model produces responses for 
the primary and one of the higher modes. See Fig. 
7 for the 5 -parameter real and imaginary compo- 
nents of 521 for the sample exhibited in Fig. 6 with 
the model for higher order modes. 

For samples 5, 6, and 7 the predicted 5 -parame- 
ter data provide realistic primary mode responses 
although the residual plots for these samples indi- 
cate that higher modes are present. With the 
higher order model described in Sec. 3.3 and the 
new solution found, the problem is resolved for the 
last three samples. The new solution specifies the 
additional parameters Pi, ... , 13^. Fig. 8 shows the 
results from the higher order model for sample 5. 
Additional work in this area suggests that this 
model is only an approximation. 




10 11 

frequency (GHz) 



10 11 

frequency (GHz) 



Fig. 3. Yttrium iron garnet predicted and observed: (a) Sf [Sn], (b) 5 [Sit]- 




10 11 

frequency (GHz) 



10 11 

frequency (GHz) 



Fig. 4. Yttrium iron garnet residual plots: (a) SR [^12], (b) 5 [Sit]- 



^ A "residual plot" denotes a plot of the difference between the predicted and observed 5-parameter data over the range of observa- 
tions. 
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Figure 5. (a) Residual plots of 5 [522] for cross-linked polystyrene (b) 3 [■^22] for an empty waveguide. 
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Figure 6. Barium titanate mix 1 (a) predicted and observed 91 [52i] (b) residual plot. 
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Figure 7. Barium titanate mix 1 (a) predicted and observed SR [52i] (b) predicted and observed 3 [Szil- 
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Figure 8. Nickel ferrite (a) predicted and observed SR [Siz] (b) predicted and observed 3 [Sn]- 



4.5 ANA Systematic Uncertainty 

The difference of the predicted 5 -parameter and 
the observed values for samples with low dielectric 
values revealed systematic uncertainty. (For the 
high dielectric samples, other sources of uncertainty 
including higher mode responses dominate the 
ANA-induced uncertainty.) Additional tests re- 
vealed that the source of the uncertainty is not re- 
lated to the material tested in the waveguide. In 
fact, uncertainties in the 5 -parameter data for the 
cross-linked polystyrene sample closely resemble 
the 5 -parameter data for an empty waveguide. At- 
tempts to further identify and remove the system- 
atic uncertainty are underway. Figure 5 contains the 
Sn data for an empty waveguide and also the resid- 
ual plot for the cross-linked polystyrene sample. 

4.6 Permittivity and Permeability Estimates 

The estimates of permeability and permittivity 
for the various samples are defined by the parame- 
ters^i, . . .A4, Bi, . . . B4 as a function of frequency. 
Several of the samples are magnetic (see Table 1). 
For these samples the permeability was allowed to 
vary, with ODRPACK required to determine the 
(1,0) value. Slight variations in the value of perme- 
ability are apparent although the deviation from the 
true value is small. The solution procedure also pro- 
vides the standard deviation and confidence inter- 
vals for each of the estimated parameters. 

4.7 Robustness of the ODRPACK Procedure 

The robustness of a mathematical procedure is its 
ability to find a locally optimal solution from a 
varied initial condition. The robustness of the entire 
permeability and permittivity procedure depends 
on the robustness of the ODRPACK procedure 
and, more significantly, the robustness of the math- 



ematical model. The existence of alternative optima 
in the mathematical model can limit the range of 
the initial conditions to produce a particular solu- 
tion. In addition, singularities in the absence of al- 
ternative optima may force the solution procedure 
to fail to determine directions of improvement and 
cause a premature termination. 

For the samples in this study, variability in the ro- 
bustness of the procedure depended on the sample 
studied. For the materials with small dielectric con- 
stant, the procedure readily determined the correct 
solution for a variety of initial conditions. For mate- 
rials with higher dielectric constant, the procedure 
often converged quickly, although the existence of 
alternative optima in the mathematical model often 
required the repeated solution with varied initial 
conditions before an acceptable solution was found. 

In particular, for the cross-linked polystyrene 
sample, initial values for fiR and €r were set to (1,0) 
and (1,0) and were constant over the entire fre- 
quency range. The solution (1,0), (2.517,0.0018) was 
found after 50 iterations of the solution procedure. 
Initial conditions of (1,0), and (4.0,0.01) failed to 
determine a solution, while altering the imaginary 
part of e^ to 1 X lO"'' resulted in the desired final 
solution once again. For material with high dielec- 
tric constants, the range was smaller. For example 
for the barium titanate mix 1, initial conditions of 
(1.0,0.0) and (200.0,1.0) produced the converged 
value of (1,0) and (269.0,1.70) in 59 iterations. 

One should note that the solution procedure can 
change the value of the estimated parameters by 
large amounts in the early iterations. Hence it is 
possible that for cross-linked polystyrene an initial 
solution farther from the desired solution may in 
fact find a solution merely because of the solution 
path taken by the procedure. 
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Figure 9. Selected plots of e^ (w) for various samples. 
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5. Conclusions and Future Directions 

The nonlinear optimization procedure using two- 
port scattering parameters determined permittivity 
and permeability for a large number of samples 
(see Fig. 9). The added capability that permits vari- 
ations in certain input parameters provides a mech- 
anism to adjust for measurement uncertainties. 
The extension of the model to include higher order 
mode responses is quite useful for high dielectric 
constant materials. 
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